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Abstract 

Coalescence of binary supermassive black holes (SBHs) would con- 
stitute the strongest sources of gravitational waves to be observed by 
LISA. While the formation of binary SBHs during galaxy mergers is al- 
most inevitable, coalescence requires that the separation between binary 
components first drop by a few orders of magnitude, due presumably to 
interaction of the binary with stars and gas in a galactic nucleus. This 
article reviews the observational evidence for binary SBHs and discusses 
how they would evolve. No completely convincing case of a bound, binary 
SBH has yet been found, although a handful of systems (e.g. interacting 
galaxies; remnants of galaxy mergers) are now believed to contain two 
SBHs at projected separations of < 1 kpc. iV-body studies of binary evo- 
lution in gas-free galaxies have reached large enough particle numbers to 
reproduce the slow, "diffusive" refilling of the binary's loss cone that is 
believed to characterize binary evolution in real galactic nuclei. While 
some of the results of these simulations - e.g. the binary hardening rate 
and eccentricity evolution - are strongly JV-dependent, others - e.g. the 
"damage" inflicted by the binary on the nucleus - are not. Luminous 
early-type galaxies often exhibit depleted cores with masses of ~ 1 — 2 
times the mass of their nuclear SBHs, consistent with the predictions of 
the binary model. Studies of the interaction of massive binaries with gas 
are still in their infancy, although much progress is expected in the near 
future. Binary coalescence has a large influence on the spins of SBHs, 
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even for mass ratios as extreme as 10 : 1, and evidence of spin-flips may 
have been observed. 

1 Introduction 

With an ever-increasing number of secure detections, supermassive black holes 
(SBHs) have evolved, in the span of a few years, from exotic possibilities to well- 
established components of galaxies. While it was understood since the 1960's 
that the energy sources of quasars must be gravitational [185], it was thirty years 
before the existence of SBHs was firmly established, through measurements of 
the Kcplcrian rise in the rotation velocity of stars or gas at the very centers 
of galactic nuclei [106]. It is now generally accepted that the formation and 
evolution of galaxies and SBHs are tightly intertwined, from the early phases 
of proto-galactic formation [200], through hierarchical build-up in CDM-likc 
cosmogonies [83], to recent galaxy mergers [150]. 

SBHs appear to be linked in fundamental ways to the dynamics of the stellar 
component in galaxies, both on large and small scales. An astonishingly tight 
correlation exists between SBH mass and the central velocity dispersion of the 
stellar component, M, <~ <r a , a w 4.5 [51]; the correlation with the velocity 
dispersion averaged over kiloparsec scales is weaker but still impressive [65, 49]. 
Similar correlations exist between SBH mass and bulge luminosity [131, 129] 
and central concentration of the light [78, 38], indicating that SBHs "know" 
about the depth of the gravitational potential well in which they live. These 
tight correlations probably reflect a degree of feedback in the growth of SBHs 
[200]. 

On small scales, SBHs are embedded in stellar cusps, parsec-scale regions 
where the stellar density increases approximately as a power law with distance 
from the SBH into the smallest resolvable radii [29, 52, 138, 66]. Faint galaxies 
have steep nuclear density profiles, p ~ r~~< ', 1.5 < 7 < 2.5, while bright galaxies 
typically have weaker cusps, 7 < 1. Steep cusps form naturally as the growth 
of the SBH pulls in stars [163]. In small dense galaxies where the star-star 
relaxation time is shorter than 10 10 yr, steep cusps may also form via collisional 
relaxation [8, 171]. Weak cusps may be remnants of strong cusps that were 
destroyed by binary SBHs during galaxy mergers; in fact the structure and 
kinematics of galactic nuclei are now believed to be fossil relics of the merger 
process [135]. 

Larger galaxies grow through the agglomeration of smaller galaxies and pro- 
togalactic fragments. If more than one of the fragments contained a SBH, the 
SBHs will form a bound system in the merger product [18, 187]. This scenario 
has received considerable attention because the ultimate coalescence of such a 
pair would generate an observable outburst of gravitational waves [209]. The 
evolution of a binary SBH can be divided into three phases [18]: 1. As the 
galaxies merge, the SBHs sink toward the center of the new galaxy via dynam- 
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ical friction where they form a binary. 2. The binary continues to decay via 
gravitational slingshot interactions [191] in which stars on orbits intersecting the 
binary are ejected at velocities comparable to the binary's orbital velocity, while 
the binary's binding energy increases. 3. If the binary's separation decreases to 
the point where the emission of gravitational waves becomes efficient at carrying 
away the last remaining angular momentum, the SBHs coalesce rapidly. 

The transition from (2) to (3) is understood to be the bottleneck of a SBH 
binary's path to coalescence, since the binary will quickly eject all stars on in- 
tersecting orbits, thus cutting off the supply of stars. This is called the "final 
parsec problem" [152]. But there are other possible ways of continuing to ex- 
tract energy and angular momentum from a binary SBH, including accretion of 
gas onto the binary system [4] or refilling of the loss cone via star-star encoun- 
ters [228, 151] or triaxial distortions [143]. Furthermore there is circumstantial 
evidence that efficient coalescence is the norm. The X-shaped radio sources 
[32] are probably galaxies in which SBHs have recently coalesced, causing jet 
directions to flip. The inferred production rate of the X-sources is comparable 
to the expected merger rate of bright ellipticals, suggesting that coalescence 
occurs relatively quickly following mergers [137]. If binary SBHs failed to merge 
efficiently, uncoalesced binaries would be present in many bright ellipticals, re- 
sulting in 3- or 4-body slingshot ejections when subsequent mergers brought in 
additional SBHs. This would produce off-center SBHs, which seem to be rare or 
non-existent, as well as (perhaps) too much scatter in the M — a and M — Ltuige 
relations [83]. 

While the final approach to coalescence of binary SBHs is not well under- 
stood, much of their dynamical effect on the surrounding nucleus takes place 
very soon after the binary forms. The binary quickly (in less than a galactic 
crossing time) ejects from the nucleus a mass in stars of order its own mass 
[177, 150] significantly lowering the central density on parsec scales. There is 
reasonably quantitative agreement between this model and the observed struc- 
ture of nuclei: the "mass deficit" - the stellar mass that is "missing" from the 
centers of galaxies, assuming that they once had steep cusps like those observed 
at the centers of faint ellipticals - is of order the black hole mass [153, 181, 77]. 

While the binary SBH model is compelling, there is still not much hard 
evidence in its support. Obscrvationally, no bona fide binary SBH (i.e. grav- 
itationally bound pair of SBHs) has definitely been detected, although there 
is circumstantial evidence (precessing radio jets; periodic outburst activity) for 
SBH binaries in a number of active galaxies, as reviewed briefly below (see 
[104] for a more complete review of this topic). But the binary SBH model has 
one great advantage: the postulated effects are accessible to observation, since 
they extend to scales of 1 — 100 pc, the distance out to which a binary SBH 
can significantly influence stellar motions. Much of the recent theoretical work 
in this field has been directed toward understanding the influence of a binary 
SBH on its stellar surroundings and looking for evidence of that influence in the 
distribution of light at the centers of galaxies. 
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Following the definition of terms and time scales in §2, we present a brief 
overview of the observational evidence for binary SBHs in §3. Interaction of 
a binary SBH with stars is discussed in §4. The possibility of multiple SBHs 
in galactic nuclei, and the implications for coalescence, are discussed in §5. §6 
summarizes TV-body work on the evolution of binary SBHs, with an emphasis 
on the question of binary wandering. Observational evidence for the destruction 
of nuclear density cusps is reviewed in §7. In some galaxies, the predominant 
source of torques leading to decay of the binary may be gas; this topic is reviewd 
in §8. Finally, the influence of binary coalescence on SBH spins is summarized 
in §9. 

2 Preliminaries 

We write mi and 771,2 for the masses of the two components of a binary SBH, 
with m 2 < mi, Q = m2/mi, and 77712 = mi + m 2 . (We also sometimes write 
M for the mass of the single SBH that forms via coalescence of two SBHs of 
combined mass m\2-) The semi-major axis of the binary's Keplerian orbit is a 
and e is the orbital eccentricity. The binary's binding energy is 

_ 6*777 1 7772 _ G^im 12 

1 1 ~ 2a ~ 2a U 
with fi — 77717772/77712 the reduced mass. The orbital period is 

/ 3 \ V 2 / \ -V 2 / \ 3/ 2 

P = ^ 7?—) = ^6 x 10 3 yr ( ( . (2) 

\Gm 12 j vio 8 aw \1 pej v ' 

The relative velocity of the two SBHs, assuming a circular orbit, is 

,_^ 658km fay . (3) 

A binary is "hard" when its binding energy per unit mass, \E\/mi2 = G/i/2a, 
exceeds <~ a 2 , where a is the ID velocity dispersion of the stars in the nucleus. 
The precise meaning of "hard" is debatable when talking about a binary whose 
components are much more massive than the surrounding stars [87, 177]. For 
concreteness, we adopt the following definition for the semi-major axis of a hard 
binary: 

„, OlS 0,, 7pc(1+9r .(^)(_^_ I )- 2 . (4) 

At distances r 3> a, stars respond to the binary as if it were a single SBH 
of mass M. The gravitational influence radius of a single SBH is defined as the 
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distance within which the force on a test mass is dominated by the SBH, rather 
than by the stars. A standard definition for r ln & is 

Thus ri n fl = 4(M/n)a,h- For an equal-mass binary, r; n fj w 16a/j, and for a more 
typical mass ratio of q = 0.1, r; n fl ss 50a^. An alternative, and often more 
useful, definition for r; n ft is the radius at which the enclosed mass in stars is 
twice the black hole mass: 

M»(r<r inf i) = 2M. (6) 

This definition is appropriate in nuclei where a is a strong function of radius; it 
is equivalent to equation (5) when the density of stars satisfies p(r) — a 2 /2irGr 2 , 
the "singular isothermal sphere," and when a is measured well outside of r; n fl. 

If the binary's semi-major axis is small enough that its subsequent evolu- 
tion is dominated by emission of gravitational radiation, then a oc — a~ 3 and 
coalescence takes place in a time t gr , where [165] 



t gr 



256F(e) G 3 fxm 2 2 ' 



This can be written 

5 Gfi 3 c 5 ( a 



tgr " 16 4 F(e) a 8 m 2 2 \a h 

3 



3.07 x 10 8 yr q 3 f 



mi2 



F{e) {l + qf \10 s MqJ ^200 km s -1 / \lQ- 2 a h 
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This relation can be simplified by making use of the tight empirical correlation 
between SBH mass and a, the "M — a relation." Of the two forms of the 
M — a relation in the literature [51, 65], the more relevant one [51] is based 
on the velocity dispersion measured in an aperture centered on the SBH, which 
is approximately the same quantity a defined above; the alternative form [65] 
defines a as a mean value along a slit that extends over the entire half-light 
radius of the galaxy. In terms of the central a, the best current estimate of the 
M — a relation is [50] 

\10 8 M Q ) y ' V200kms-V V ' 
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with a = 4.86 ± 0.43. Combining equations (8) and (9) and setting M = m\i, 
F = 1 gives 



3 / \ -3.14 / \ 4 

Coalescence in a Hubble time (~ 10 10 yr) requires a < 0.05a/j for an equal-mass 
binary and a < 0.15a^ for a binary with q = 0.1. Inducing a SBH to decay from 
a separation a w ah ~ 10° pc to a separation such that £ gr < 10 10 yr is called 
the "final parsec problem" [152]. Much of the theoretical work on massive black 
hole binary evolution has focussed on this problem. 

3 Observations of Binary Supermassive Black 
Holes 

3.1 External Galaxies 

If a binary SBH is defined as two SBHs separated by a distance a < ah, then 
no completely convincing example of such a binary has yet been found. Here 
we briefly review the small set of cases in which clear evidence is seen for two, 
widely separated SBHs in a single system ("dual SBHs"), as well as the still- 
circumstantial evidence for true binary SBHs. For a more complete review of 
this topic, see [104]. 

3.1.1 Dual SBHs 

Figure 1 shows what was probably the first clear example of two SBHs in one 
"system," in this case a pair of interacting galaxies near the center of the galaxy 
cluster Abell 400. The associated radio source 3C75 consists of a pair of twin 
radio lobes originating from the radio cores of the two galaxies; the projected 
separation of the cores is ~ 7 kpc [161]. Such double-jet systems are expected 
to be rare given the small fraction of giant elliptical galaxies that are associated 
with luminous radio sources. 

"Binary" quasars are common but most are believed to be chance projections 
or lensed images [155, 102]. Among the binary quasars for which lensing can 
be ruled out, the smallest projected separation belongs to LBQS 0103-2753 at 
z = 0.85, with an apparent spacing between centers of 2.3 kpc [95]. However 
the two quasar spectra show a Az of 0.024 suggesting a chance projection. 

Galaxies in the late stages of a merger are the most plausible sites for dual 
SBHs and many of these exhibit double nuclei in the optical or infrared [80, 
23] . However few show unambiguous evidence of AGN activity in both nuclei, 
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Figure 1: 20 cm VLA image of the radio source 3C 75 in the cluster of galaxies 
Abell 400. The image consists of two, twin-jet radio sources associated with each 
of two elliptical galaxies. The jets bend and appear to be interacting. The pro- 
jected separation of the radio cores is about 7 kpc. Image courtesy of NRAO/AUI 
and F. N. Owen et al. 

indicative of SBHs. One clear exception is NGC 6240 (Fig. 2), for which both 
nuclei exhibit the flat X-ray spectra characteristic of AGNs [105]. The projected 
separation is 1.4 kpc. Another likely case is Arp 299 [12]. 

Interestingly, there are no known dual SBHs with separations below ~ 1 
kpc, even though a 1 kpc separation would be resolvable to distances of several 
hundred Mpc. 

3.1.2 Evidence for Binary SBHs 

Many active galaxies exhibit periodic variability with periods of days or years, 
consistent with the orbital periods of true binary SBHs having a < a^. Un- 
doubtedly the clearest example is OJ 287, a "blazar," i.e. an active galaxy in 
which the jet is believed to be orientated nearly parallel to the line of sight, at 
z = 0.306. Optical variability of OJ 287 has been recorded since 1890 [175, 206] 
and has a strict period of 11.86 yr (~ 9 yr in the galaxy's rest frame); the last 
major outburst was observed (on schedule) in 1994. The outbursts are gener- 
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Figure 2: Chandra X-ray image of the starburst galaxy NGC 624-0, showing the 
two nuclear sources. Projected separation of the nuclei is about 1.4 kpc. Image 
courtesty of NASA/CXC/MPE/S.Komossa et al. 

ally double-peaked with the peaks separated by about a year; the second peak 
is accompanied by enhanced radio emission. Models to explain the periodicity 
usually invoke a second SBH with q < 0.1. In one class of model, the variability 
reflects true changes in the source luminosity due to variations in the accre- 
tion rate as the smaller SBH passes through the accretion disk surrounding the 
larger SBH [201, 117, 213]. In these models, the observed variability period 
is equal to the binary orbital period, and the binary orbit is highly eccentric 
(e sa 0.7), implying a relatively short (< 10 5 yr) time scale for orbital decay via 
gravitational radiation. The lag between primary and secondary peaks may be 
due to the time required for the disturbance induced by the passage through 
the accretion disk to propagate down the jet [213]. Alternatively, the luminosity 
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variations may reflect changes in the jet direction resulting from precession of 
the accretion disk, the latter induced by torques from the second SBH [96]. In 
this model, the binary orbital period is much less than 9 yr, and the secondary 
maxima could be due to a "nodding" motion of the accretion disk [96]. 

Many other examples of variability in AGN at optical, radio and even TeV 
energies are documented [223], with periods as short as ~ 25 days [85]. Indeed 
evidence for variability has even been claimed for the Milky Way SBH, at radio 
wavelengths; the ostensible period is 106 days [231]. However none of these 
examples exhibits as clear a periodicity as OJ287. 

Table 3.1.2 gives a list of active galaxies for which periodic variability has 
been claimed. 



Source 


Period (yr) 


Reference 


Mkn 421 


23.1 


[121] 


PKS 0735+178 


14.2 


[46] 


BL Lac 


14 


[47] 


ON 231 


13.6 


[122] 


OJ 287 


11.9 


[175] 


PKS 1510-089 


0.92 


[225] 


Sgr A* 


0.290 


[231] 


3C 345 


10.1 


[230] 


AO 0235+16 


5.7 


[180] 


3C 66A 


0.175 


[110] 


Mkn 501 


0.065 


[85] 


3C 273 


0.0026 


[224] 



Table 1: Sources with periodic variation in the nuclear emission 

Radio lobes in active galaxies provide a fossil record of the orientation his- 
tory of the jets powering the lobes. Many examples of sinusoidally or helically 
distorted jets are known, and these observations are often interpreted via a bi- 
nary SBH model. The wiggles may be due to physical displacements of the 
SBH emitting the jet (e.g., [188]) or to precession of the larger SBH induced by 
orbital motion of the smaller SBH (e.g., [186]). In the radio galaxy 3C 66B, the 
position of the radio core shows well-defined elliptical motions with a period of 
just 1.05 yr [204], implying t gr < 10 3 yr. 

About a dozen radio galaxies exhibit abrupt changes in the orientation of 
their radio lobes, producing a "winged" or A-shaped morphology [115]. While 
originally interpreted via a precession model [36], a more likely explanation is 
that the SBH producing the jet has undergone a spin flip, due perhaps to capture 
of a second SBH [137, 233]. 

A number of quasars show the peaks of their broad emission lines at very 
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different redshifts from their narrow emission lines, or two displaced emission 
line peaks, which might be attributed to orbital motion of the SBHs associated 
with the line emitting regions [63, 64, 203, 232]. This interpretation has fallen 
out of favor however since the candidate systems do not show the predicted 
radial velocity variations [37]. 

A number of other possibilities exist for detecting binary SBHs, including 
the use of space interferometers to measure the astrometric reflex motion of 
AGN photocenters due to orbital motion of the jet-producing SBHs [219]; mea- 
surement of periodic shifts in pulsar arrival times due to passage of gravitational 
waves from binary SBHs [123]; and, of course, direct detection of gravitational 
waves by space-based interferometers. 

3.2 Limits on the Binarity of The Milky Way Black Hole 

The likely longevity of binary SBHs motivates the question whether the clos- 
est and best-studied SBH at the center of the Milky Way galaxy is a binary. 
Monitoring of the proper motion of stars orbiting the SBH has led to a precise 
measurement of its mass, M~4x 10 6 M© [69, 193]. The Milky Way SBH is 
coincident with the compact (< 1 AU) radio source Sagittarius (Sgr) A*. If the 
Milky Way SBH were a binary, the radio source would probably be associated 
with the more massive of the two binary components. Limits on the masses 
of the components could be placed by measuring astrometric reflex motion of 
the radio source relative to distant quasars [6, 184]. Such measurements have 
recovered the magnitude of the solar reflex motion in the galaxy but have so far 
yielded no evidence for a binary SBH. The most recent upper limits on the mass 
of a binary companion of Sgr A* are M 2 < 10 4 M Q for binaries with semimajor 
axes 10 3 AU < a < 10 5 AU [183]. This places any companion that may exist 
in the class of "intermediate-mass" black holes (IBHs). The parameter space of 
SBH-IBH binaries at the Galactic center is illustrated in Figure 3. 

IBHs have been suggested as a possible explanation for ultraluminous X-ray 
sources; however their existence is not widely accepted. It has been suggested 
that the center of the Milky Way is a place where IBHs might naturally form 
via the runaway merging of massive stars in the young, dense star clusters ( [59] 
and references therein). Two such clusters, the Arches and the Quintuplet, are 
presently located in the Galactic center region. The segregation of massive stars 
to the cluster center accelerates the "core collapse" in which the stellar density at 
the center of the cluster increases drastically. Collapse time can be shorter than 
the life time of the most massive stars; in this case runaway stellar coalescence 
ensues resulting in the formation of a supermassive star at the cluster center. 
If the star survives mass loss through winds and avoids exploding as a pair- 
instability supernova, it collapses to form an IBH [220]. Dynamical friction in 
the background stellar cusp of the Galactic bulge subsequently drags the IBH 
toward the SBH until two black holes form a hard binary. This process might 
explain the puzzling presence of early-type stars [67, 69] deep inside the sphere 
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Figure 3: A crude illustration of the parameter space for a SBH-IBH binary at 
the Galactic center. Assuming a circular orbit around a SBH of 3 x 1O 6 M0 ; a 
IBH with mass Mibh and semi-major axis a can be ruled out by measurement 
of an astrometric wobble of the radio image of Sgr A* . The shaded regions show 
the detection thresholds for astrometric resolutions o/O.Ol, 0.1 and 1 milliarc- 
seconds, respectively, assuming a monitoring period of 10 years. The dashed 
lines indicate coalescence due to gravitational radiation in 10 6 and 10 7 years, 
respectively (From [84], see also [229]). 

of influence of the SBH at the Galactic center ([84], but see [99]). 
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4 Interaction of Binary Black Holes with Stars 



4.1 Dynamics of a Massive Binary in a Fixed Stellar Back- 
ground 

Stars passing within a distance ~ 3a of the center of mass of a hard binary 
undergo a complex interaction with the two black holes, followed almost always 
by ejection at velocity ~ \J /V m i2^bin, the "gravitational slingshot" [191]. Each 
ejected star carries away energy and angular momentum, causing the semi- 
major axis, eccentricity, orientation, and center-of-mass velocity of the binary 
to change and the local density of stars to drop. If the stellar distribution is 
assumed fixed far from the binary and if the contribution to the potential from 
the stars is ignored, the rate at which these changes occur can be computed 
by carrying out scattering experiments of massless stars against a binary whose 
orbital elements remain fixed during each interaction [89, 187, 87, 88, 13, 148, 
177, 133, 134]. Figure 4 shows an example of field star velocity changes in a set 
of scattering experiments. 

Consider an encounter of a single field star of mass rn* with the binary. Long 
before the encounter, the field star has velocity v with respect to the center 
of mass of the field star-binary system and its impact parameter is p. A long 
time after the encounter, the velocity v of the field star attains a constant value. 
Conservation of linear momentum implies that the change 6V in the velocity of 
the binary's center of mass is given by 

SV = -H^-5v. (11) 

TO12 

The velocity change results in a random walk of the binary's center-of-mass 
momentum, as discussed in more detail below. The energy of the field star- 
binary system, expressed in terms of pre-encounter quantities, is 

_ 1 o 1 T .o Gm\mi 

E = -m*t; + -ro 12 V --^— 

2 V m 12 J 2a K ' 

with V = — (m Jr /mi 2 )v the initial velocity of the binary's center of mass and 
do the binary's initial semi-major axis. After the encounter, 

1 f ^ m+ \ 2 Gmim 2 



E = -m* 1 + — \v A (13) 

2 \ mi2 / 2a 



and E = E , so that 



m±(v 2 



Gra\m2 \ m\2 
m±(v 2 - vp 
Gm\m<2, 



(14) 



12 



i 1 1 1 1 1 1 1 1 1 1 1 1 r 




Figure 4: Distribution of field star velocity changes for a set of scattering ex- 
periments in which the field star's velocity at infinity relative to the binary was 
v/o = O.bVbin&x- The binary's mass ratio was 1 : 1, and the orientation of the 
binary's orbital plane with respect to the x-axis was varied randomly between 
the scattering experiments. Each plot represents 5 x 10 4 scattering experiments 
within some range of impact parameters [pi,P2] in units of a. (a) [6,10] (b) 
[2,4] (c) [0.6,1] (d) [0.4,0.6]. Solid lines in (a) and (b) are the distributions 
corresponding to scattering off a point-mass perturber. In (c) and (d), the mean 
of this distribution (which is very narrow) is indicated by the arrows. The grav- 
itational slingshot is apparent in the rightward shift of the Sv values when p is 
small, due to the randomization of ejection angles (from [133]). 

Averaged over a distribution of field-star velocities and directions, equation 
(4.1) gives the binary hardening rate (d/dt){l/a) [87, 88, 148, 177]. 

The angular momentum of the field star-binary system about its center of 
mass, expressed in terms of pre-encounter quantities, is 

C a = mjl + —) £ + M 4o (15) 

where £o = P v o and £bo = £bo 1 1 1 with Cb the binary's orbital angular momen- 
tum. Conservation of angular momentum during the encounter gives 

M V m 12/ 
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(16) 



Changes in \£{,\ correspond to changes in the binary's orbital eccentricity e via 
the relation e 2 = 1 — £ 2 /Gm\2H 2 a [148, 177]. Changes in the direction of if, 
correspond to changes in the orientation of the binary [134]. 

The results of the scattering experiments can be summarized via a set of 
dimensionless coefficients H, J, K, L, ... which define the mean rates of change 
of the parameters characterizing the binary and the stellar background. These 
coefficients are functions of the binary mass ratio, eccentricity and hardness but 
are typically independent of a in the limit that the binary is very hard. The 
hardening rate of the binary is given by 

d ^ (17) 



dt \a 

with p and a the density and ID velocity dispersion of stars at infinity. The 
mass ejection rate is 

dWfa) = Jmi2 (18) 

with M e j the mass in stars that escape the binary. The rate of change of the 
binary's orbital eccentricity is 

de K. (19) 



dln(l/a) 

The diffusion coefficient describing changes in the binary's orientation is 

{M , )=L m^Gpa 
mi2 a 

with the stellar mass. Brownian motion of the binary's center of mass is 
determined by the coefficients A and C which characterize the Chandrasekhar 
diffusion coefficients at low V: 

(At;,,) = -AV, 

(Av^) = 1 -(Avl) = C. (21) 

The mean square velocity of the binary's center of mass is (V 2 ) = C/2A. 

The binary hardening coefficient H reaches a constant value of <~ 16 in the 
limit a <C Ofc, with a weak dependence on q [87, 148, 177]. In a fixed background, 
equation (17) therefore implies that a hard binary hardens at a constant rate: 

1 1 Go 

— w H (t — th) , t>t h , a{t h ) = a h . (22) 

a(t) a h a 
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This is sometimes taken as the definition of a "hard binary." The time to reach 
zero separation is 



t-t h = 



a 1 2cr 3 



2g(l + g) 2 



.3 



HGp a h HG 2 pp 



H 



G 2 pmi2 



5.2 x 10 5 yr q(l + qf 
( mi2 



200 km s 




1O 3 M pc- 3 



P 



) 



-l 



X 




(23) 



Orbital shrinkage would occur quite rapidly in the environment of a galactic 
nucleus if the properties of the stellar background remained fixed. 

However if the binary manages to shrink to a separation at which t gr < 10 10 
yr, the changes it induces in its stellar surroundings will be considerable. The 
mass ejected by the binary in decaying from ah to a gr is given by the integral 
of equation (18): 



Figure 5 shows M e j as a function of the mass ratio q for a = 200 km s -1 and 
various values of t gr . The mass ejected in reaching coalescence is of order m\2 for 
equal- mass binaries, and several times 777,2 when 771,2 <C m\. A SBH that grew to 
its current size through a succession of mergers should therefore have displaced 
a few times its own mass in stars. If this mass came mostly from stars that 
were originally in the nucleus, the density within ri n fl would drop drastically 
and the hardening would stop. Without some way of replenishing the supply of 
stars (and in the absence of other mechanisms for extracting angular momentum 
from the binary, e.g., torques from gas clouds; cf. §8), decay would stall at a 
separation much greater than a gr . 

In the case of extreme binary mass ratios, 7772 <C mi, the assumption that 
all stars passing a distance ~ a from the binary will be ejected is likely to be 
incorrect; many such stars will pass through the binary system without being 
appreciably perturbed by the smaller black hole. The concept of "ejection" 
in extreme binary mass ratios may be misleading; since the typical amount of 
energy transferred in stellar pericenter passages is small, most of the slingshot 
stars are not ejected from the nucleus and remain bound to the larger black 
hole. Figure 5 should thus be interpreted with caution in the regime 7772 <C 777,1. 
The extreme-mass-ratio regime is poorly understood but deserves more study 
in view of the possibility that intermediate-mass black holes may exist having 
masses of - 10 3 M Q [43]. 

Changes in the binary's orbital eccentricity (eq. 19) are potentially impor- 
tant because the gravity wave coalescence time drops rapidly as e — > 1 (cq. 
8). For a hard binary, scattering experiments give K(e) w K e(l — e 2 ), with 
K ks 0.5 for an equal-mass binary [148, 177]. The dependence of K on m^/mi 
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Figure 5: Mass ejected by a decaying binary, in units of mu = mi + m 2 (solid 
lines) or m 2 (dashed lines), calculated by an integration of equation (24), with 
the coefficient J (a) taken from [177]. Curves show mass that must be ejected in 
order for the binary to reach a separation where the emission of gravitational 
radiation causes coalescence on a time scale of 10 10 yr (lower), 10 9 yr (middle) 
and 10 8 yr (upper). 

is not well understood and is an important topic for further study. The implied 
changes in e as a binary decays from a = ah to zero are modest, Ae < 0.2, for 
all initial eccentricities. 

4.2 Evolution in an Evolving Background 

The scattering experiments summarized above treat the binary's environment 
as fixed and homogeneous. In reality, the binary is embedded at the center of an 
inhomogeneous and evolving galaxy, and the supply of stars that can interact 
with it is limited. 

In a fixed spherical galaxy, stars can interact with the binary only if their 
pericenters lie within <~ IZxa, where 1Z is of order unity. Let Li c = [E — Q(lZa)] w 

V2Gmi27^a, the angular momentum of a star with pericenter TZa. The "loss 
cone" is the region in phase space defined by L < Li c . The mass of stars in the 
loss cone is 
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. j dE J ^ dL 2 4n 2 f(E,L 2 )P(E,L 2 ) 

w &ir 2 Gm 12 m*lla J dEf(E)P rad (E). (25) 

Here P is the orbital period, / is the number density of stars in phase space, and 
N(E, L 2 )dEdL is the number of stars in the integral-space volume defined by dE 
and dL. In the final line, / is assumed isotropic and P has been approximated 
by the period of a radial orbit of energy E. An upper limit to the mass that is 
available to interact with the binary is <~ M; c (a^), the mass within the loss cone 
when the binary first becomes hard; this is an upper limit since some stars that 
are initially within the loss cone will "fall out" as the binary shrinks. Assuming 
a singular isothermal sphere for the stellar distribution, p oc r~ 2 , and taking the 
lower limit of the energy integral to be &(a h ), equation (25) implies 

M lc (a h ) w 3Kp. (26) 

We can compute the change in a that would result if the binary interacted with 
this entire mass, by using the fact the mean energy change of a star interacting 
with a hard binary is <~ 3Gp/2a [177]. Equating the energy carried away by 
stars with the change in the binary's binding energy gives 



2 a 2 
or 

, /a h \ AM 9Kp _ q 
V a J mi2 mi2 (1 + q) 

if AM is equated with Mi c . Only for very low mass ratios (q < 1CP 3 ) is this 
decay factor large enough to give t gr < 10 10 yr (eq. 8), but the time required 
for such a small black hole to reach the nucleus is likely to exceed a Hubble time 
[132]. Hence even under the most favorable assumptions, the binary would not 
be able to interact with enough mass to reach gravity-wave coalescence. 

But the situation is even worse than this, since not all of the mass in the loss 
cone will find its way into the binary. The time scale for the binary to shrink is 
comparable with stellar orbital periods, and some of the stars with r per i s=s 
will only reach the binary after a has fallen below ~ a h . We can account for 
the changing size of the loss cone by writing 

dM f x 1 dM lc 



-dE 

dt J Eo (t)P(E) dE 



= 8ir 2 Gm 12 m*lla(t) / fi(E)dE, (29) 
JE (t) 

where M{t) is the mass in stars interacting with the binary and fi(E) is the 
initial distribution function; setting P(E ) = t reflects the fact that stars on 
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orbits with periods less than t have already interacted with the binary and been 
ejected. Combining equations (27) and (29), 

- ] w 247r 2 "KGm* / U{E)dE. (30) 

a J JE {t) 

Solutions to equation (30) show that a binary in a singular isothermal sphere 
galaxy stalls at ah/a w 2.5 for m2 = mi, compared with a/j/a w 10 if the full 
loss cone were depleted (equation 28). In galaxies with shallower central cusps, 
decay of the binary would stall at even greater separations. 




4.3 Collisional Loss-Cone Replenishment 

A binary black hole depletes its loss cone very quickly since stars within the loss 
cone need only a few close encounters with the binary to be ejected. Whether 
the binary can continue to exchange energy with stars depends on the efficiency 
with which stars are re-supplied to the loss cone. 

The most commonly invoked mechanism for loss cone re-filling is two-body 
scattering of stars. A small angular momentum perturbation, for instance from 
a passing star, can deflect a star with L > Li c into the loss cone. This process 
has been studied in detail in the context of scattering of stars into the tidal 
disruption sphere of a single black hole [56, 119, 27]. The basic equations are 
similar in the case of scattering into a binary SBH, except that the critical an- 
gular momentum increases by a factor <~ ya/rt , where r t is the tidal disruption 
radius. (Other differences are discussed below.) If the binary parameters are 
assumed fixed, a steady-state flow of stars into the loss cone will be achieved 
on roughly a two-body relaxation time scale Tr, and the distribution function 
near Li c will have the form 

where R is a scaled angular momentum variable, R = L 2 /Ll(E), L C (E) is the 
angular momentum of a circular orbit of energy E, and / is the distribution 
function far from the loss cone, assumed to be isotropic. The mass flow into the 
central object is m* / T{E)dE, where 

T(E)dE = 4n 2 L 2 c (E) U — lim } —I—dE. (32) 

V ' cV ; \7 v r R^o 2R J ln(l/i& c ) 1 ; 

The quantity in brackets is the orbit-averaged diffusion coefficient in R. 
A crude estimate of the collisional re-supply rate is given by 
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[56], where M*(r) is the mass in stars within radius r and r cr it is the critical 
radius at which stars scatter into the loss cone in a single orbital period; beyond 
r cr it, the diffusion rate drops rapidly with radius. Estimates based on simple 
galaxy models give r cr ;t w 10 — 100a. To get an idea of the scattering rate, we 
consider the nucleus of the Milky Way. Equation (4) gives ah ~ 0.32(1 + q) _1 
pc. Assuming a — ah, q = 0.1 and r cr ; t w 30(1^ gives r cr ; t w 1.1 pc. The 
mass within this radius is ~ 3 x 1O 6 M [67] and the relaxation time at this 
radius, assuming stars of a solar mass, is <~ 2 x 10 9 yr. The scattered mass 
over 10 10 yr is then <~ 10 7 M Q . This is comparable to the mass of the Milky 
Way SMBH, M « 3 x 10 6 M© [193], but the scattering rate would drop as the 
binary shrinks, suggesting that scattered stars would contribute only modestly 
to refilling of the loss cone. In more massive galaxies, the nuclear density is 
lower, relaxation times are longer, and collisional refilling would be even less 
important. A more detailed calculation of the collisional refilling rates in real 
galaxies [228] concludes that few if any binary SBHs could reach coalescence via 
this mechanism. The author took the presently observed luminosity profiles of 
galaxies as initial conditions. The stellar density profiles must have been steeper 
before the binary SBH formed [150], leading to substantially more rapid decay 
early in the life of the binary. 

Another criticism of standard loss cone theory is its assumption of a quasi- 
steady-state distribution of stars in phase space near Li c [151]. This assumption 
is appropriate at the center of a globular cluster, where relaxation times arc 
much shorter than the age of the universe, but is less appropriate for a galactic 
nucleus, where relaxation times almost always greatly exceed a Hubble time 
[44]. (The exceptions are the nuclei of small dense systems like the bulge of 
the Milky Way.) The distribution function f(E,L) immediately following the 
formation of a hard binary is approximately a step function, 



much steeper than the ~ In L dependence in a collisonally relaxed nucleus (equa- 
tion 31). Since the transport rate in phase space is proportional to the gradient 
of / with respect to L, steep gradients imply an enhanced flux into the loss cone. 
Figure 6 shows the evolution of N(E, R) at a single E assuming that the loss 
cone is empty initially within some R[ c = Li c /L C (E) and that N(E,R,t = 0) 
is a constant function of R outside of Ri c . (The loss cone boundary is as- 
sumed static; in reality it would shrink with the binary.) Also shown is the 
collisionally-relaxed solution of equation (31). The phase-space gradients decay 
rapidly at first and then more gradually as they approach the steady-state solu- 
tion. The total mass consumed by the binary, shown in the lower panel of Fig. 
6, is substantially greater than would be computed from the steady-state the- 
ory, implying greater cusp destruction and more rapid decay of the binary. This 
time-dependent loss cone refilling might be particularly effective in a nucleus 
that that continues to experience mergers or accretion events, in such a way 




L > Lie 
L < Lie, 



(34) 
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Figure 6: (a) Slices of the density N(E,R,t) at one, arbitrary E, recorded, 
from left to right, at 10°, 10 1 , 10 2 7 10 3 , and 10 4 Myr (solid curve). Initially, 
N(E,R,t) — for R < R ic and N(E, R, t) = const for R > R lc . We also 
show the equilibrium solution of equation (31) (dot-dashed curve), (b) The total 
number of stars consumed by the loss cone as a function of time (solid curve). 
The scale has been set to galaxy M32 with initial separation between the MBHs 
of 0.1 pc. (From [151]) 

There are other differences between loss cones around single and binary 
SBHs. A star that interacts with a massive binary generally remains inside the 
galaxy and is available for further interactions. In principle, a single star can 
interact many times with the binary before being ejected from the galaxy or 
falling outside the loss cone; each interaction takes additional energy from the 
binary and hastens its decay. Consider a simple model in which a group of N 
stars in a spherical galaxy interact with the binary and receive a mean energy 
increment of (AE). Let the original energy of the stars be E . Averaged over a 
single orbital period P(E), the binary hardens at a rate 

d (Gm 1 m 2 \ N(AE) 



dt 



V 2a ) ~ m * P(E) • 



m*— — — — . (35) 
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In subsequent interactions, the number of stars that remain inside the loss cone 
scales as Lf c cx a while the ejection energy scales as ~ a -1 . Hence N(AE) cx 
a 1 a _1 cx a . Assuming the singular isothermal sphere potential for the galaxy, 
one finds 

t m*N(AE) t-t h 



a h 1 . M , 
— w 1 + In 

a(t) mi2 



(36) 



2 M a 2 P(£ ). 

[151]. Hence the binary's binding energy increases as the logarithm of the time, 
even after all the stars in the loss cone have interacted at least once with the 
binary. Re-ejection would occur differently in nonspherical galaxies where angu- 
lar momentum is not conserved and ejected stars could miss the binary on their 
second passage. However there will generally exist a subset of orbits defined by 
a maximum pericenter distance < a and stars scattered onto such orbits can 
continue to interact with the binary. 

As these arguments suggest, the long-term evolution of a binary SBH due 
to interactions with stars may be very different in different environments. (We 
stress that the presence of gas may substantially alter this picture; cf. §8.) 
There are three characteristic regimes (we stress that the presence of gas may 
substantially alter this picture; cf. section 8) [151]. 

1. Collisional. The relaxation time Tr is shorter than the lifetime of the 
system and the phase-space gradients at the edge of the loss cone are given 
by steady-state solutions to the Fokker-Planck equation. The densest galactic 
nuclei may be in this regime. Resupply of the loss cone takes place on the time 
scale associated with scattering of stars onto eccentric orbits. The decay time 
of a binary SBH scales as |a/d| ~ m^ 1 <~ N with N the number of stars. In 
the densest galactic nuclei, collisional loss cone refilling may just be able to 
drive a binary SBH to coalescence in a Hubble time. For sufficiently small Tr, 
scattering refills the loss cone in less than an orbital period ("full loss cone") 
and the decay follow a -1 ~ t. A-body studies are typically in this regime, as 
discussed below. 

2. Collisionless. The relaxation time is longer than the system lifetime and 
gravitational encounters between stars can be ignored. The low-density nuclei 
of bright elliptical galaxies are in this regime. The binary SBH quickly interacts 
with stars whose pericenters lie within its sphere of influence; in a low-density 
(spherical or axisymmetric) nucleus, the associated mass is less than that of the 
binary and the decay tends to stall at a separation too large for gravitational 
wave emission to be effective. However evolution can continue due to re-ejection 
of stars that lie within the binary's loss cone but have not yet escaped from the 
system. In the spherical geometry, re-ejection implies \a/a\ ~ (1 + t/t )/a, 
leading to a logarithmic dependence of binary hardness on time. Re-ejection 
in galactic nuclei may contribute a factor of ~ a few to the change in a over a 
Hubble time. 

3. Intermediate. The relaxation time is of order the age of the system or 
somewhat longer. While gravitational encounters contribute to the re-population 



21 



Tabic 2: Physical Regimes for Long-Term Decay of Massive Black Hole Binaries 



Form of Decay 


Regime 


a -1 oc const 


Collisionless 


cx t/N 


Collisional (diffusion) 


cx t + const 


Collisional (full loss cone) 


cx ln(l + t/to) + const 


Rc-ejcction 



of the loss cone, not enough time has elapsed for the phase space distribution 
to have reached a collisional steady state. Most galactic nuclei are probably in 
this regime. The flux of stars into the loss cone can be substantially higher than 
predicted by the steady-state theory, due to strong gradients in the phase space 
density near the loss cone boundary produced when the binary SBH initially 
formed. This transitory enhancement would be most important in a nucleus 
that continues to experience mergers or infall, in such a way that the loss cone 
repeatedly returns to an unrelaxed state with its associated steep gradients. 

Table 1 summarizes the different regimes. The evolution of a real binary SBH 
may reflect a combination of these and other mechanisms, such as interaction 
with gas. There is a close parallel between the final parsec problem and the 
problem of quasar fueling: both requre that of order 10 8 M Q be supplied to the 
inner parsec of a galaxy in a time shorter than the age of the universe. Nature 
clearly accomplishes this in the case of quasars, probably through gas flows 
driven by torques from stellar bars. The same inflow of gas could contribute 
to the decay of a binary SBH in a number of ways: by leading to the renewed 
formation of stars which subsequently interact with the binary; by inducing 
torques which extract angular momentum from the binary; through accretion, 
increasing the masses of one or both of the SBHs and reducing their separation; 
etc. 

4.4 Non-axisymmetric nuclei 

The estimates made above were based on spherical models of nuclei. The total 
number of stars in a full loss cone can be much larger if the nucleus is flattened 
and axisymmetric [124], when only one component of the angular momentum is 
conserved. In very flattened nuclei (with ellipticities e <~ 0.5), single emptying 
of an initially full loss cone can in some cases be sufficient to drive the binary to 
coalescence [228]. However loss cone dynamics can be qualitatively different in 
non-axisymmetric (triaxial or bar-like) potentials, since a much greater number 
of stars may be on "centrophilic" - box or chaotic - orbits which take them 
arbitrarily near to the SBH(s) [159, 68, 196, 212, 228]. Stars on centrophilic 
orbits of energy E experience pericenter passages with r por j < d at a rate ~ 
A(E)d [143]. If the fraction of stars on such orbits is appreciable, the supply 
of stars into the binary's loss cone will remain essentially constant, even in 
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the absence of collisional loss-cone refilling. Such models need to be taken 
seriously given recent demonstrations [169, 90, 170] that galaxies can remain 
stably triaxial even when composed largely of centrophilic orbits. Furthermore 
imaging of galaxy centers on parsec scales reveals a wealth of features in the 
stellar distribution that are not consistent with axisymmetry, including bars, 
nuclear spirals, and other misaligned features [221, 164, 40]. 

The total rate at which stars pass within a distance IZa of the massive binary 

is 

^ w KaJ A(E)M c (E)dE (37) 

where M c (E)dE is the mass on centrophilic orbits in the energy range E to 
E + dE. In a nucleus with p ~ r~ 2 , the implied feeding rate into a radius r; n fl 
is roughly 

M* « Tc^r (38) 

w 2500M Q yr^T ( ») , (39) 

°' y 1 V200 kms 3 y V ' 

where f c is the fraction of stars on centrophilic orbits. If this rate were main- 
tained, the binary would interact with its own mass in stars in a time of only 
<~ 10 5 yr, similar to the decay time estimated above (eq. (23) for a binary in 
a fixed background. In fact, the feeding rate would decline with time as the 
centrophilic orbits were depleted. Solving the coupled set of equations for a(t) 
and M c (t), one finds that at late times, the binary separation in a p cx r~ 2 
nucleus varies as [170] 

Comparison with Table 1 shows that this is the same time dependence as for 
the "full loss cone" regime of spherical nuclei. Placing just a few percent of a 
galaxy's mass on centrophilic orbits is sufficient to overcome the final parsec 
problem and induce coalesence, if the stellar density profile is steep and if the 
chaotic orbits are present at all energies. This example is highly idealized, but 
shows that departures from axial symmetry in galactic nuclei can greatly affect 
the rate of decay of a binary SBH. 



5 Multiple Black Hole Systems 

If binary decay stalls, an uncoalesced binary may be present in a nucleus when 
a third SBH, or a second binary, is deposited there following a subsequent 
merger. The multiple SBH system that forms will engage in its own gravitational 
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slingshot interactions, eventually ejecting one or more of the SBHs from the 
nucleus and possibly from the galaxy and transferring energy to the stellar 
fluid. 

If the infalling SBH is less massive than either of the components of the pre- 
existing binary, m 3 < (mi,m 2 ), the ultimate outcome is likely to be ejection of 
the smaller SBH and recoil of the binary, with the binary eventually returning to 
the galaxy center. The lighter SBH is ejected with a velocity roughly 1 /3 the rel- 
ative orbital velocity of the binary [191, 92], and the binary recoils with a speed 
that is lower by m-sjimx + m.2). Each close interaction of the smaller SBH with 
the binary increases the latter's binding energy by (E/E) w Arris/ (mi + m<i) 
[89]. If m 3 > mi or m 3 > m 2 , there will most often be an exchange interaction, 
with the lightest SBH ejected and the two most massive SBHs forming a binary; 
further interactions then proceed as in the case m.3 < (mi, m.2). 

During the three-body interactions, both the semi-major axis and eccentric- 
ity of the dominant binary change stochastically. Since the rate of gravity wave 
emission is a strong function of both parameters (E oc a~ 4 (l — e 2 )~ 7 / 2 ), the 
timescale for coalescence can be enormously shortened. This may be the most 
promising way to coalesce SBH binaries in the low-density nuclei of massive 
galaxies, where stalling of the dominant binary is likely. 

This process has been extensively modelled using the PN2.5 approximation 
to represent gravitational wave losses [166] and assuming a fixed potential for 
the galaxy [214, 147, 215]. In these studies, there was no attempt to follow 
the pre-merger evolution of the galaxies or the interaction of the binary SBHs 
with stars. In two short non-technical contributions (submissions for the IEEE 
Gordon Bell prizes in 2001 and 2002), J. Makino and collaborators mention 
two iV-body simulations of triple SBH systems at the centers of galaxies using 
the GRAPE-6, and (apparently) a modified version of NBODY1. Relativistic 
energy losses were neglected and the SBH particles all had the same mass. 
Plots of the time evolution of the orbital parameters of the dominant binary 
show strong and chaotic eccentricity evolution, with values as high as 0.997 
reached for short periods. Such a binary would lose energy by gravity wave 
emission very rapidly, by a factor ~ 10 8 at the time of peak e compared with a 
circular-orbit binary with the same semi-major axis. 

In a wide, hierarchical triple, m 3 <C (mi, m.2), the eccentricity of the domi- 
nant binary oscillates through a maximum value of <~ yjl — 5 cos 2 i/3, | cosi| < 



y/3/5, with i the mutual inclination angle [107]. One study [21] estimates that 

the coalescence time of the dominant binary in hierarchical triples can be re- 
duced by factors of <~ 10 via the Kozai mechanism. 

If the binary SBH is hard when the third SBH falls in, the ejected SBH can 
gain enough velocity to escape the galaxy. If the three masses are comparable, 
even the binary can be kicked up to escape velocity. One study [217] estimates 
(based on a very simplified model of the interactions) that the recoil velocity of 
the smallest SBH is larger than galactic escape velocities in 99% of encounters 
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and that the binary escapes in 8% of encounters. Thus a signficant fraction of 
nuclei could be left with no SBH, with an offset SBH, or with a SBH whose 
mass is lower than expected based on the M — a or M — Lbuige relations. 

There is a need for simulations of multiple-SBH systems that include both 
gravitational loss terms, accurate (regularized) interactions between the SBH 
particles, and the interactions of SBH particles with stars. 

6 TV-Body Studies of Binary Black Hole Evolu- 
tion 

The interaction of a massive binary with point perturbers at the center of a 
galaxy is a straightforward problem for TV-body simulation. In principle, TV- 
body studies can reveal both the long-term evolution of the binary, as well as 
the effect of the binary on its stellar surroundings. The latter can be compared 
with observed nuclear density profiles as a test of the theory (§7). 

Unfortunately, unless great care is taken, TV-body studies are likely to give 
misleading results. This follows from the result (Section 4.3) that time scales 
for two-body scattering of stars into the binary's loss cone are of order 10 10 yr 
or somewhat longer in real galaxies. In TV-body simulations, relaxation times 
are shorter by factors of <~ TV/10 11 than in real galaxies, hence the long-term 
evolution of the binary is likely to be dominated by spurious loss cone refilling, 
wandering of the binary, and other noise-driven effects. 

TV-body studies are most useful at characterizing the early stages of binary 
formation and decay, or simulating the disruptive tidal effects of a SBH on the 
nucleus of an infalling galaxy. Due to algorithmic limitations - primarily the 
difficulty of integrating galaxy models with high central concentrations - most 
such studies [35, 76, 126, 125, 158, 157, 86, 25, 128] have been based on galaxy 
models with unrealistically large cores. 

Figure 7 is from the first [136] TV-body simulation of galaxy mergers in 
which the pre-merger galaxies contained power-law nuclear cusps as well as 
massive particles representing the SBHs. These simulations were run using 
GADGET [202], a tree code with inter-particle softening, and were not able to 
accurately follow the formation and decay of the massive binary. The SBH in 
the larger galaxy was found to tidally disrupt the steep cusp in the infalling 
galaxy, producing a remnant with only slightly higher central density than that 
of the giant galaxy initially. This result helps to explain the absence of dense 
cusps in bright galaxies [55] , and suggests that the central structure of galaxies 
can only be understood by taking into account the destructive influence of SBHs 
on the stellar distribution during mergers. Additional results, using a similar TV- 
body code and a variety of mass ratios for the merging galaxies, were reported 
in [141]. 

Had these simulations been extended to longer times using a more accurate 
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Figure 7: final density profiles from a set of 10:1 merger simulations in which 
each galaxy contained a black hole (a-d) and in which neither galaxy contained 
a black hole (e-g) [136]. The four thin curves in each frame correspond to four 
different pre-merger orbits, (a), (e) Space density of stars initially associated 
with the secondary galaxy; thick curves are the initial density profile, (b), (f) 
Space density of stars initially associated with the primary galaxy; thick curves 
are the initial density profile, (c), (g) Space density of all stars. Lower thick 
curves are the initial density profile of the primary galaxy, and upper thick 
curves are the superposition of the initial density profiles of the primary and 
secondary galaxies. Lines of logarithmic slope — 1 and —2 are also shown, (d), 
(h) Logarithmic slope of the surface density profiles of the merger remnants. 
Thick curves correspond to the initial primary galaxy. 

N-hoAy code, the massive binary would have ejected stars via the gravitational 
slingshot and lowered the central density still more. This was first demonstrated 
[150] in an iV-body study that used a tree code for the early stages of the 
merger, and NB0DY6, a high-precision, direct-summation code [1], for the later 
stages, when the binary separation fell below the tree code's softening length. 
The pre-merger galaxies had steep, p ~ r~ 2 density cusps and the mass ratio 
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was 1:1. These simulations were continued until the binary separation had 
decayed by a factor of ~ 10 below an- The initially steep nuclear cusps were 
converted to shallower, p ~ profiles shortly after the SBH particles had 
formed a hard binary; thereafter the nuclear profile evolved slowly toward even 
shallower slopes as the massive binary ejected stars. As Fig. 5 illustrates, the 
stellar density around the binary drops very quickly after the binary reaches a 
separation a s=s a^. 
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Figure 8: Lagrangian radii around each of the two SBH particles in an equal- 
mass merger simulation [150]. $From bottom to top, the radii enclose 10~ 4 , 
10~ 3 - 5 , 10~ 3 , 10~ 2 - 5 , 10~ 2 7 1CT 1 - 5 and 1CT 1 in units of the mass of one galaxy 
before the merger. The binary becomes "hard" at t w 11, and very rapidly heats 
the surrounding stellar fluid, lowering the local density. 

The hardening rate of the binary in these simulations was found not to be 
strongly dependent on the number of particles. This result was subsequently 
shown [151] to be due to the small N: stars were resupplied to the loss cone 
via collisions at a higher rate than they were being kicked out by the binary, 
ensuring a continuous supply of stars and allowing the binary to continue to 
shrink. While a qualitatively similar evolution may take place in some galaxies 
- for instance, loss cones in non-axisymmetric potentials can be continuously 
repopulated by stars on centrophilic orbits (§4.4) - collisional loss cone refilling 
is very unlikely to achieve anything like a full loss cone except in very small, 
dense galaxies. The long-term evolution of the binary in almost all published 
TV-body simulations are therefore not representative of what one would expect 
in real galaxies. 
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Figure 9: Evolution of the binary semi-major axis (a) and hardening rate (b) in 
a set of high accuracy N-body simulations; the initial galaxy model was a low- 
central- density Plummer sphere ([20]). Units are G — M ga i — 1, E = —1/4, 
with E the total energy, (a) Dashed lines are simulations with binary mass Mi = 
Mi = 0.005 and solid lines are for Mi = Mn = 0.02, in units where the total 
galaxy mass is one. (b) Filled(open) circles are for Mi — M% = 0.005(0.02). 
Crosses indicate the hardening rate predicted by a simple model in which the 
supply of stars to the binary is limited by the rate at which they can be scattered 
into the binary's influence sphere by gravitational encounters. The simulations 
with largest (Mi, M 2 ) exhibit the nearly A -1 dependence expected in the "empty 
loss cone" regime that is characteristic of real galaxies. 



q(E) = (41) 



The relevant dimensionless parameter is 

_ (5 J? 
% 

where 5 J is the change over one radial period in the angular momentum of 
a star on a low-J orbit, and J; c is the angular momentum of an orbit at the 
edge of the binary's loss cone. A value q(E) 3> 1 implies that the loss cone 
orbits at energy E are re-populated at a much higher rate than they are de- 
populated by the binary, and the loss cone remains nearly full. A value g« 1 
implies that the loss cone is essentially empty, and repopulation must talke 
place diffusively, as stars scatter in from J > J; c . In real galaxies, A is large 
and S J is small, implying q <C 1. Achieving q -C 1 in A-body simulations 
requires large particle numbers, and/or a model for the galaxy that has an 
unrealistically low central density, so that the star-star relaxation time is long. 
Figure 9 shows a recent set of A-body simulations that does both [20]. The 
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massive binary was embedded in a Plummcr [168] galaxy model, which has 
a core radius comparable to its half-mass radius; this model is very different 
from real galaxies but its very low degree of central concentration implies a 
long relaxation time and low rate of collisional loss-cone refilling. Large particle 
numbers were achieved, without sacrificing accuracy, by running the simulation 
on a parallel GRAPE cluster. The TV-dependence of the binary's hardening 
rate is clear; in the simulations with binary mass Mi = M 2 = 0.02M ffa ;, the TV- 
dependence of the hardening rate is s = (d/dt)(l/a) w A~ 0,8 , almost as steep 
as the A -1 dependence predicted for an "empty" loss cone [151]. A similar 
study [128], based on King-model galaxies, found a similar result. 

The Plummer-model initial conditions used in the simulations of Figure 9 
were identical to those adopted in two other A-body studies based on a more 
approximate A-body code [178, 25]. Contrary to Figure 9, Chatterjee, Hcrn- 
quist & Loeb (2003) found that the binary hardening rate "saturated" at values 
of A > 2 x 10 , remaining constant up to A w 4 x 10 5 . They speculated that 
this was due to a kind of Brownian-motion-mediated feedback, in which the bi- 
nary maintains a constant supply rate by modulating the local density of stars. 
However no supporting evidence for this model was presented; for instance, it 
was not demonstrated that the central density was actually regulated by the 
binary, or that the amplitude of the Brownian wandering increased with A in 
the manner postulated. Furthermore these authors provided no plots showing 
the claimed A-dependence of the hardening rate. Chatterjee et al.'s conclusion, 
that "a substantial fraction of all massive binaries in galaxies can coalesce within 
a Hubble time," is not substantiated by the more accurate A-body simulations 
shown in Figure 9. 

While Brownian motion probably does affect the decay rate of binaries in A- 
body simulations [151], it is doubtful that the effect is significant in real galaxies. 
The Brownian velocity of single black holes is found in A-body integrations to 
be [114] 

l -M{V 2 ) * (42) 

where a 1 is the ID, mean square stellar velocity within a region r < 0.5r^ 
around the black hole (and includes the influence of the black hole on the stellar 
motions), and is the stellar mass. In the case of the Milky Way black hole, 
equation (42) implies V lms w 0.2 km s _1 (assuming to* = M Q ) and an rms 
displacement of < 0.1 pc. Brownian motion of a massive binary is larger than 
that of a single black hole, but only by a modest factor [133, 150, 128]. The rms 
displacement of a binary from its otherwise central location would therefore be 
very small in a real galaxy, probably even less than the separation between the 
two components of the binary. 

The goal of A-body studies is to simulate binary evolution in galaxies 
with realistic density profiles, and with large enough particle numbers that re- 
population of the binary's loss cone takes place diffusively, as in real galaxies. 
Two avenues are open for making further progress in this area. 1. One could 
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Figure 10: Results from a set of N-body integrations of a massive binary in a 
galaxy with a p ~ 7-~°- D density cusp [205]. Each curve is the average of a set of 
integrations starting from different random realizations of the same initial condi- 
tions, (a) Evolution of the "mass deficit" (Eq. J i 3), i.e. the mass in stars ejected 
by the binary. For a given value of binary separation a, the mass deficit is nearly 
independent of particle number N , implying that one can draw conclusions from 
observed mass deficits about the binary that produced them, (b) Evolution of bi- 
nary eccentricity. The eccentricity evolution is strongly N -dependent and tends 
to decrease with increasing N , suggesting that the eccentricity evolution in real 
binaries would be modest. 

combine a Monte-Carlo treatment of stellar encounters [57, 58] with a lookup 
table, derived from scattering experiments, of energy and angular momentum 
changes experienced during close passages of stars to the binary [178, 133]. Such 
a hybrid algorithm would allow one to adjust the degree of collisionality at will 
and record the effects on both the binary's evolution, and the influence of the 
binary on the stellar distribution. This approach would be difficult to generalize 
to non-spherical geometries however. 2. A straightforward iV-body approach is 
also feasible, but particle numbers in excess of ~ 10 7 are required [151]. Such 
large particle numbers are just now becoming feasible for direct-summation N- 
body codes, by combining the GRAPE accelerator boards [127] with a parallel 
architecture [33]. Indeed just such an approach was used for the integrations of 
Figure 9, although the galaxy models in that study were rather unphysical. 

Figure 10 shows a promising early step in this direction. The initial galaxy 
models had p ~ r ~ - 5 density cusps; integrations were carried out on a GRAPE- 
6 computer, limiting the total particle number to ~ 256K, but the motion of 
the black hole binary (of mass Mi = M2 = 0.005M go ;) and nearby stars was 
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carried out using the Mikkola-Aarseth chain rcgularization algorithm [145, 146, 
3]. Because of the models' higher central density and limited particle numbers, 
the binary's loss cone was only partially empty, q > 1, and the TV-dependence of 
the hardening rate was shallower than expected in real galaxies, (d/dt)(l/a) ~ 
jy-o.4 pig urc io (a) shows that the "damage" inflicted by the binary on the 
nucleus is not strongly dependent on N, as expected (cf. Eq. 18). This is an 
encouraging result since it implies that one can hope to learn something definite 
about pre-existing binaries by comparing TV-body simulations with observations 
of the centers of current-day galaxies (§7). On the other hand, Figure 10 (b) 
suggests that the evolution of the binary's eccentricity is strongly TV-dependent. 
This may explain the rather disparate results on eccentricity evolution in past 
TV-body studies [150, 86, 2]. 

Much progress on this problem is expected in the next few years. 

7 Evidence for Cusp Destruction 

A potentially powerful constraint on models of binary SBH evolution is the 
observed central density structure of galaxies. Figure 5 shows that a massive 
binary must eject of order its own mass in reaching a separation at which t gr < 
10 10 yr if m2 ~ mi, or several times 777,2 if 777,2 "C mi. These numbers should be 
interpreted with caution since: (1) binaries might not decay this far - they may 
stall - or the final stages of decay might be driven by gas dynamics rather than 
energy exchange with stars; (2) the definition of "ejection" used in Figure 5 
is escape of a star from an isolated binary, and does not take into account 
the confining effect of the nuclear potential; (3) the effect of repeated mergers 
on nuclear density profiles, particularly mergers involving very unequal-mass 
binaries, is poorly understood. Neverthless, even the initial formation of a hard 
binary displaces a mass of order m 2 (Fig. 8). 

The luminosity profile data can probably be used to rule out one model of 
binary evolution. In a "collisionless" galaxy (Table 1), the binary's loss cone 
never refills, and decay of the binary would stall. The binary carves out a "hole" 
in both phase space and configuration space; the radius of the latter would 
be ~ 3ah [233]. While central minima may have been seen in the luminosity 
profiles of a few galaxies [113], these are likely due to dust obscuration, and the 
great majority of galaxies show a clearly rising stellar density into radii < rj n fl. 
The non-existence of true "cores" suggests either that some degree of loss-cone 
refilling occurs, or that the final decay of the binary takes place via a more 
efficient process than ejection of stars. 

Nevertheless there is a well-defined trend for the central densities of bright 
galaxies to decrease with increasing luminosity [52, 138, 44, 66, 79]. Nuclear 
densities in elliptical galaxies and spiral bulges with My < — 20 follow p ~ r~ 7 , 
7 < 1, while in fainter spheroids, 1 < 7 < 2.5. A natural interpretation is that 
the brightest galaxies - which presumably formed via one or more mergers - 
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have experienced more cusp destruction than fainter galaxies. (An alternative 
possibility, discussed below, is that the nuclei in faint galaxies re-formed after 
being destroyed.) 

In practice, this hypothesis is difficult to test since it requires knowledge of 
the pre-merger density profiles. A reasonable guess is that all galaxies originally 
had steep power-law density cusps, since these are generic in the faintest galaxies 
known to harbor SBHs. For instance, both M32 and the bulge of the Milky Way 
have p ~ r -1 ' 5 at r < r infl and p ~ r~ 2 just outside [112, 67]. 
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Figure 11: Observed surface brightness profile of NGC 3348. The dashed line is 
the best-fitting Sersic model to the large-radius data. Solid line is the fit of an 
alternative model, the "core- Sersic" model, which fits both the inner and outer 
data well. The mass deficit is illustrated by the area designated "depleted zone" 
and the corresponding mass is roughly 3 x 10 8 M Q [77]. 

The "mass deficit" [153] is defined as the difference in integrated mass be- 
tween the observed density profile and the primordial (pre-merger) profile. For 
instance, if the primoridal profile is a power law of index 70 inward of some 
radius rj,, then 
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Mass deficits in samples of bright elliptical galaxies were computed in three 
recent studies [153, 181, 77]. In the first two studies, the authors assumed power- 
laws of various slopes for the pre-merger profiles, and found (M^ef/M,) w 1 for 
7o = 1.5 with M. the current mass of the SBH. The latter study made use of the 
fact that the light profiles of bright galaxies show an abrupt downward deviation 
relative to a Scrsic [197] profile fit to the outer regions (Fig. 11). Mass deficits 
inferred in this study were slightly larger, M de f/M. w 2. 1 These numbers are 
within the range predicted by the binary SBH model, particularly given the 
uncertainties associated with the effects of multiple mergers. 

In small dense galaxies, a destroyed cusp would be expected to re-form via 
the Bahcall-Wolf [8, 171] process, on a timescale of order the star-star relaxation 
time measured at r; n fl. This time is of order 10 9 yr in the Milky Way bulge and 
the nucleus of M32. This may be the explanation for the steep power-law profiles 
observed at the centers of these galaxies. Alternatively, the steep cusps may be 
due to star formation that occurred after the most recent merger [97]. 

More rigorous tests of the binary SBH model will require a better under- 
standing of the expected effect of massive binaries on stellar density profiles. 
As discussed above, while the best current TV-body simulations suggest p <~ r _1 
following binary formation [150], the simulations are dominated by noise over 
the long term. 

A number of other processes could compete with binary SBHs in the de- 
struction of nuclear density cusps. A population of three or more SBHs in a 
galactic nucleus would undergo a complicated set of close encounters resulting 
ultimately in coalescence and/or ejection of some or all of the SBHs (§5). In the 
process, the stellar background would be heated and a mass of order five times 
the combined mass in SBHs removed [139]. This model reproduces the observed 
time dependence of core radii in globular clusters very well [142] but its rele- 
vance to galactic nuclei is less clear; the model requires binary coalescence times 
long enough that an uncoalesccd binary is present when a third SBH falls in 
[217]. If a binary SBH does eventually coalesce, the gravitational radiation car- 
ries a linear momentum leading to a recoil of the coalesced hole [19, 53]. Recoil 
velocities are estimated to be as large as ~ 400 km s _1 [48, 140], although with 
considerable uncertainty. A SBH ejected from a galactic nucleus with a velocity 
of ~ 10 2 km s _1 would quickly fall back to the center, but its displacement and 
infall would heat the stellar fluid and lower its density. Figure 12 shows the 
effects of ejection on nuclear density profiles. Mass deficits produced by this 
mechanism can be comparable in amplitude to those predicted by the binary 
SBH model. 

A major focus of future work should be to calculate the evolution of p(r) as 
predicted by the various scenarios for binary decay discussed in this article. 

lr The author of this study presents his mass deficits as significantly smaller than those 
found in the earlier studies. However he bases his comparison on values of M^cf computed 
exclusively using 70 = 2. 
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Figure 12: Effect on the nuclear density profile of SBH ejection. The initial 
galaxy model (black line) has a p ~ r _1 density cusp, (a) Impulsive removal 
of the SBH. Tick marks show the radius of the black hole 's sphere of influence 
Hnfl before ejection. A core forms with radius ~ 2n n fl. f&j Ejection at velocities 
less than escape velocity. The black hole has mass 0.3% that of the galaxy; the 
galaxy is initially spherical and the black hole 's orbit remains nearly radial as it 
decays via dynamical friction. The arrow in this panel marks n n fl in the initial 
galaxy. [140]. 

8 Interaction of Binary Black Holes with Gas 

Interstellar gas might play an important role in the dynamical evolution of 
a binary SBH. Interactions with gas complement interactions with the stellar 
environment (§4) and with other SBHs (§5). Any gas situated close to a binary is 
disturbed by the SBHs and exerts gravitational force on them, thereby affecting 
their orbit. Furthermore, if SBH coalescence is accompanied by the presence of 
gas, an observable electromagnetic afterglow might follow coalescence. 

The collisional, dissipative nature of interstellar gas gives rise to a behavior 
fundamentally different from that of the point-mass dynamics of stellar systems. 
It is natural to distinguish between two classes of flows in dynamical systems 
containing gas. In hot flows the gas temperature is comparable to the virial 
temperature of the system, while in cold flows the gas temperature is signifi- 
cantly below the virial temperature. The virial temperature can be defined as 
T v ir = GMi2fim p /2ak 7 where fi is the mean particle mass in units of the proton 
mass m p , and k is the Boltzmann constant. The prototype of a hot flow is the 
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spherical, "Bondi" accretion onto a single black hole, in which the accreting gas 
is supported by pressure against free infall toward the accretor. The prototype 
of a cold flow is a thin disk, in which the gas is rotationally supported against 
infall. Even in hot flows rotational support is realized close to the accretor when 
the gas has nonzero net angular momentum (e.g. [109]). 

The angular momentum barrier is central to SBH formation theories. Any 
model for how material is channeled into an accreting black hole must describe 
the mechanism by which angular momentum is removed from the material. 
Whatever this mechanism may be, it is expected that it operates universally 
during the epoch in which SBHs grew to their present masses by rapidly ac- 
creting material onto pre-existing black hole "seeds." This is also the period 
when galaxy merging peaks [97, 83, 222, 216]. While still elusive to astronomi- 
cal probes due to severe obscuration [190], the nuclei of merging galaxies, which 
are also the sites for the formation of binary SBHs [18], are expected to contain 
the largest concentration of dense gas anywhere in the universe. The inevitable 
abundance of gas motivates an inquiry into the role of gas dynamics as an al- 
ternative to stellar dynamics in the process of SBH coalescence. Some of the 
mechanisms that remove angular momentum of interstellar gas and thus chan- 
nel it into the neighborhood of SBHs include the torquing of gas flow by the 
rapidly-fluctuating potential of merging galaxies [144] and by nested stellar bars 
[199], angular momentum transport by hydrodynamical turbulence that might 
be driven by the onset of self-gravity [198, 61, 74] or by supernovae embed- 
ded within a large-scale toroidal circumnuclear flow [218], angular momentum 
extraction by magnetohydrodynamical turbulence [11] or by magnetic braking 
[22], and more speculatively, by Rossby vortex instabilities [118]. 

Astronomical observations offer abundant evidence for both hot and cold 
gas flows in the immediate vicinity of SBH candidates. The origin and the 
dynamical impact of the two classes of gas flow are distinct and are discussed 
here separately. 

8.1 Interaction with Hot Gas 

Hot gas permeates the interstellar space in galaxies and the intergalactic space 
in groups of galaxies and galaxy clusters. Virial temperatures range between 
10 6 K — 10 8 K and the hot gas is almost completely ionized. Primordial and 
secondary sources contribute to the pool of hot has. During the early stages 
of galaxy formation, intergalactic space contains partially ionized gas inherited 
from the pregalactic, early universe. Hydrogen recombines at redshifts z <~ 1000 
and is reionized at redshifts z <~ 10 by the radiation emitted by the earliest 
structures. The partially ionized gas cools within the confining gravitational 
potential of dark matter halos and filaments. Cold gas accelerates toward the 
halos' centers of gravity and is shock-heated to about the virial temperature. 
Some of the coldest inflowing gas escapes heating by accreting along narrow 
channels that reach deep inside the primary halo. Cooling times in the halo 
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centers where the gas is the densest are short compared to the dynamical time 
and thus most of the primordial gas is consumed in starbursts on a dynamical 
time scale. 

Tenuous gas that remains after the cooling time has exceeded the dynamical 
time in the nascent galaxy might still be plentiful enough to feed a massive black 
hole growing at an Eddington-limited rate. The residual number density at the 
radius of influence of the SBH is 
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where T is the virial temperature of the galaxy, \i is the average atomic mass in 
units of the proton mass, A is the cooling function [30], and we have employed 
the M - a relation (Eq. 9) to relate the virial temperature to the black hole 
mass. The thermal stability limit could in principle be exceeded if the gas kept 
at the Compton temperature by a continuum flux from an unobscured AGN 
[45]. 

This so-called "cooling flow model of quasar fueling" [26, 160] is however 
plagued by many problems ([108] and references therein). Most of the gas 
left over from star formation might be blown out by the mechanical feedback 
associated with the radiative and mechanical output of the accreting massive 
black hole [200, 100, 156]. A small amount of angular momentum in the gas 
results in circularization and settling into an accretion disk. This disk may be 
susceptible to fragmentation, thereby converting most of the gas mass into stars 
and effectively cutting off the supply of gas to the SBH [208] . 

The geometry of the flow of a hot, magnetized gas near a binary black hole is 
unknown. Assuming spherical, non-rotating accretion, the time scale on which 
the hot gas is captured by the SBH is 
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where /b ~ 1 — 10 is a numerical factor that depends on the equation of state 
of the gas. 

If a binary black hole is present, gravitational torques from the gas induce 
decay of the binary's semi-major axies on approximately the same time scale. 
This crude estimate is based on an analogy with binary-star interactions: the 
binary must eject of order its own mass in stars to decay an e-folding in separa- 
tion. Hot gas torquing the binary might be ejected in an outflow and thus the 
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actual rate at which gas is accreting onto individual binary components might 
be severely suppressed compared to the accretion expected in an isolated black 
hole. 

Galactic nuclei also contain hot gas produced by secondary sources. For 
example, observations with the Chandra X-ray Observatory have revealed ten- 
uous (n ps 10—100 cm -3 ), hot (T ps 1 keV) plasma within a parsec of the 
~4x 10 6 M Q Milky Way SBH [7]. This plasma is being generated by the nu- 
merous massive, evolved stars in the galactic region [67] through stellar wind 
and supernova activity. Since its temperature is higher than the virial, most 
(> 99%) of the plasma escapes the neighborhood of the SBH [176]. While the 
hot gas densities in active galaxies might be transiently larger than that at the 
Galactic center, the tendency of the hot plasma to escape the neighborhood 
of the SBH reduces the likelihood that large quantities of virialized gas would 
remain enmeshed with the binary's orbit long enough to affect its dynamical 
evolution. 

Recently, Escala et al. [41, 42] carried out smoothed particle hydrodynamical 
(SPH) simulations of binary point masses interacting with a massive, spherical 
cloud of hot gas initially centered on the binary. Gravitational drag from the 
gas induces decay in the binary's orbit. The relevance of spherical, hot ini- 
tial conditions is contingent on the astrophysical plausibility that a compressed 
accumulation of hot gas comparable in mass to the SBH can be sustained. 

8.2 Interaction with Cold Gas 

The specific angular momentum of a cold flow might easily exceed that of the 
binary. The gas then tends to settle into rotationally supported, geometrically 
thin rings and disks (recall that "cold" gas is colder than the virial temperature 
but can be hot enough to be ionized). 

Observations offer abundant evidence for the presence of dense gas in galactic 
nuclei. Thin, Keplerian molecular disks on scales 0.1 pc — 0.5 pc have been 
seen in the water maser emission in the nuclei of Seyfert galaxies [154, 60, 
82]. The Galactic nucleus contains a 4 x 1O 6 M0 black hole surrounded by a 
~ (10 4 — 10 5 )M Q molecular gas torus at distances > 1 pc from the SBH [94]. 
Compact stellar disks on scales > 20 pc, which are fossil evidence of past gas 
circularization, are evident in the nuclei are of many galaxies [167]. Massive 
accretion disks must be present in quasars and the Narrow-Line Seyfert I nuclei 
to account for what appears to be rapid accretion onto the central SBHs in 
these systems. However, the structure of these disks at radii comparable to the 
size of a hard SBH binary is unknown. The principal challenges to constructing 
extended disk models are the instabilities related to incomplete ionization and 
the susceptibility to gravitational fragmentation [103, 198, 149, 75]. 

If a disk surrounding a binary SBH is initially inclined with respect to the 
binary's orbital plane, the quadrupole component of the binary's gravitational 
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potential causes differential precession in the disk at the rate [111] 
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which results in a warping of the disk. As in the Bardeen-Petterson mechanism 
[15, 173, 194], the warp either dissipates, or smears around the binary, resulting 
ultimately in a nearly axisymmetric disk in the binary's orbital plane. 

Interest in co-planar, circumbinary disks stems from their ability to extract 
a binary's angular momentum via a form of tidal coupling. Two interrelated 
questions might be posed: What is the response of a circumbinary disk to the 
binary's tidal forcing? and: How does such a disk affect the evolution of the 
binary's orbit? 

Existing attempts to answer these questions have employed ad hoc models 
for the form of the binary-disk torque coupling [172, 93], or have been restricted 
to binaries with components of very unequal mass [4] , where an array of neigh- 
boring Lindblad resonances facilitate binary-disk coupling [72], much like the 
coupling between a massive planet and its natal gas disk [73]. Early numer- 
ical simulations of circumbinary disks with nearly equal masses [5], however, 
suggested that the disks are truncated exterior to the resonances, which was 
interpreted as a consequence of a collisionless nonlinear parametric instability 
[189, 39] . Fluid dynamical theory of circumbinary disk truncation is still lacking. 

In a circular binary the outer Lindblad resonances (OLR) are located at 
radii r m = (1 + 1/m) 2 / 3 , where m = 1, 2, ... is the order in the decomposition of 
the binary's gravitational potential into multipoles: 

oo 

<p(r, 0)=^ ip m (r) cos[to((9 - Q hin t)}. (47) 

The outermost OLR is located at r w 1.6a. The resonances are radii in the 
disk where the natural, epicyclic frequency of radial oscillations in the disk is an 
integer multiple of the rate at which a packet of disk gas receives tidal "kicks" 
by the binary. The forcing near a resonance, as well as at a radius where 
surface density in the disk exhibits a large gradient, excites nonaxisymmetric 
propagating disturbances, or "density waves," in the disk. 

The gravitational potential of eccentric binaries contains low-frequency com- 
ponents that are absent in circular binaries. These low-frequency components 
activate resonances located at larger radii than in the circular case, and might 
lead to mutual excitation and reinforcement of the binary and the disk ec- 
centricities [162, 71]. Many extrasolar planets, which are thought to form in 
circumstellar disks, are notably eccentric 2 , suggesting that dynamical coupling 
between a binary point mass (a star and a planet, or a pair of black holes) 
and a gas disk is conducive to eccentricity growth. The observed circumbinary 

2 http: / /cxoplanets.org 
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disks in young stellar binaries such as GG Tau, which are typically eccentric, 
are truncated at radii a few times the semimajor axis [130], which lends support 
to this hypothesis. Eccentricity in SBH binaries accelerates coalescence due to 
gravitational wave emission (Eq. 7) and might be detectable in gravitational 
wave trains. 

Density waves transport angular momentum outward through the circumbi- 
nary disk. Angular momentum flux carried by the waves is extracted from the 
binary's angular momentum. The binary experiences a negative torque equal 
and opposite to the total angular momentum flux transferred to the disk. The 
location of the inner edge of the disk reflects a balance between the angular 
momentum flux deposited into the disk, and the angular momentum flux trans- 
ported through the disk by another, possibly viscous mechanism. Wave mo- 
mentum is deposited into the disk material via a form of dissipative damping. 
The location in the disk where the waves are damped can be separated by many 
wavelengths from the location where they are excited. The damping could take 
place in the nonlinear steepening and the breaking of wave crests [192, 179]. In 
marginally optically thick disks, radiation damping might also play a role [24]. 
Yet another form of damping could be due to the dissipation of wave shear if the 
disk is strongly viscous [207] . The amplitude of the density waves is a steeply 
decreasing function of the radius of excitation. The amplitude is diminished if 
the waves are nonlinear at excitation and damp in situ, but then one expects 
the inner edge to recede where in situ damping shuts off. 

The intricate and insufficiently understood nature of binary-disk interactions 
calls for grid-based hydrodynamical simulations with a shock-capturing capabil- 
ity. The necessity that the radial wavelength, which is smaller than the vertical 
scale height of the disk, be resolved by multiple cells, places severe demands 
on the computational resources, especially if a three dimensional representation 
of the disk is required. It should also be noted that the radiative and thermal 
structure of accretion disks around single SBHs are not adequately understood 
on any radial scale. 

As a binary's semimajor axis decreases due to stellar, gas dynamical, or 
gravitational radiation processes, a circumbinary disk's inner edge spreads in- 
ward viscously while maintaining constant edgc-to-scmimajor axis ratio, e.g., 
r 'edge /a ~ 2. In the final stages of the gravitational radiation-driven inspiral, 
however, the time scale on which the semimajor axis decays becomes shorter 
than the viscous time scale, and the disk can no longer keep up with the bi- 
nary, resulting in binary-disk detachment. On the relevant length scales the disk 
might be dominated by radiation pressure and the electron scattering opacity; 
the structure and the stability of such disks is an active research area [210]. 
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9 Spin Evolution during Mergers 

Coalescence of a binary black hole results in a spinning remnant. 3 Angular 
momentum conservation implies 

Si + S 2 + L orb = S + J rad (48) 

where Si and S 2 are the spin angular momenta of the two SBHs just before 
the final plunge, L or b is the orbital angular momentum of the binary before the 
plunge, S is the spin of the resulting black hole, and J ra d is the angular momen- 
tum carried away by gravitational waves during and after the coalescence [54]. 
The simplest case to treat is extreme mass ratio mergers, q = m 2 /mi <C 1, for 
which the binary can be described as a test particle of mass m 2 = m orbiting a 
black hole of mass mi = M ^> m, and both S 2 and J ra d can be ignored. The 
change in the larger hole's spin is computed by adding the smaller hole's energy 
and orbital angular momentum at the last stable orbit (LSO). The latter varies 
from £lso / m = \/l2M for circular equatorial orbits around a non-spinning hole 
to Llso = M(9M) for prograde(retrograde) orbits around a maximally-spinning 
hole, Si = M 2 . The much larger value of Llso m the case of retrograde capture 
implies that a rapidly-rotating hole will typically spin down if capture occurs 
from random directions [34, 70, 226, 227, 91]. The change in spin assuming 
q < 1 is 

5a = q(-2a + L LS o.,z) + 0(q 2 ) (49) 

where a = |Si|/M 2 , L z is the orbital angular momentum parallel to S, and 
L = L/mM. The first term in equation (49) describes conservation of spin 
angular momentum of the larger hole as its mass grows, a cx M~ 2 , while the 
second term describes the increase in spin due to torquing by the smaller body. 

The change in spin after a single coalescence is illustrated in Fig. 13 as 
a function of q and initial spin; the upper (lower) curves represent prograde 
(retrograde) captures from equatorial orbits, and the dashed lines are for capture 
over the pole. The bias toward spin-down is evident; retrograde capture from 
the equatorial plane produces a nearly (but never completely) non-spinning 
remnant when q s=a 2.5a, q < 0.23, and rapid final rotation (a > 0.9) requires 
both a large initial spin and a favorable inclination. On the other hand, if 
the larger hole is slowly rotating initially, a < 0.5, mass ratios q > 0.3 always 
result in spin-up. The oft-repeated statement that "mergers spin down black 
holes" reflects a preconception that SBHs are likely to be formed in a state of 
near-maximal rotation (e.g., [14, 62]). 4 

Successive mergers from random directions with fixed q (i.e. secondary mass 
grows proportionately to primary mass) lead to a steady-state spin distribution 

3 In this section, we set G = c = 1. 

4 Note the error in Figure 1 of Hughes & Blandford (2003), which shows the change in 
spin for mergers with mass ratio q = 0.5: the darkest contour on that plot should be labelled 
a = 0.5, not a = 0. 
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Figure 13: Final spin a of a remnant black hole in terms of its original spin, for 
mass ratios q = 0.1 (red), 0.3 (green) and 0.5 (blue). The change in spin was 
computed using the test-particle approximation for £lso [16, 91, 227]. Upper 
(lower) curves correspond to prograde (retrograde) capture from the equatorial 
plane; dashed curves are for capture over the pole. Capture of a low-mass sec- 
ondary is likely to spin down the larger hole unless the latter is slowly rotating 
initially. Capture of a massive secondary results in spinup unless infall is nearly 
retrograde or the original spin is large. 



N(a) that is uniquely determined by q. For small q, this distribution can be 
derived from the Fokker-Planck equation [91]: 

N(a)da w 7V a 2 e" 3a2/2a ™ s da, a rms « 1.58^. (50) 

Figure 14 shows N(a) for various values of q, computed via Monte-Carlo experi- 
ments (not from the Fokker-Planck equation) using the test-mass approximation 
for Llso- The Gaussian form of equation (50) is seen to be accurate only for 
q < 0.1. For q > 1/8, the distribution is skewed toward large spins. 

Accurate calculation of spin- up during a merger of comparably massive black 
holes requires a fully general-relativistic numerical treatment. Adopting vari- 
ous approximations for the radius of the innermost stable circular orbit (ISCO) 
for comparably massive binaries [28, 17, 81, 31], and assuming that mass and 
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Figure 14: Steady-state spin distributions produced by successive capture from 
random directions at fixed mass ratio q, for q = (1/32,1/16,1/8,1/4,1/2). 
Curves were generated using Monte- Carlo experiments based on the test-particle 
approximation, q <C 1; hence the curve for q = 1/2 should be viewed as illustra- 
tive only. 

angular momenta are conserved during coalescence, gives a remnant spin in 
equal-mass mergers of a w 0.8 — 0.9. Baker et al. [9, 10] present full numer- 
ical calculations of equal-mass mergers with and without initial spins. In the 
absence of initial spins, 3% of the system's mass-energy and 12% of its angular 
momentum are lost to gravitational radiation, and the final spin is a w 0.72. Co- 
alescence of initially spinning holes from circular orbits in the equatorial plane 
yields a w 0.72 + 0.32s with s the initial spin parameter of the two holes (as- 
sumed equal); Baker et al. considered initial spins in the range —0.3 < s < 0.2, 
where negative/positive values indicate spins aligned/countcraligned with the 
orbital angular momentum. Extrapolating this result toward s = 1 suggests 
that prograde mergers of black holes with initial spins s > 0.85 will result in a 
maximally- spinning remnant. 

Confronting these predictions with observation is problematic for a number 
of reasons: merger histories of observed SBHs are not known, SBH spins are 
difficult to determine observationally, and other mechanisms, such as gas accre- 
tion, can act efficiently to spin up SBHs [14]. However there is circumstantial 
evidence that mergers played a dominant role in determining the spins of at least 
some SBH. If SBH spins were the product of gas accretion, the jets in active 
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galaxies should point nearly perpendicularly to the disks of their host galaxies 
[182]. In fact, there is almost no correlation between jet direction and galaxy 
major axis in Seyfert galaxies [211, 101, 195]. Among the possible explanations 
[101] for the misalignment, perhaps the most natural is that SBH spins were 
determined by the same merger events that formed the bulge, long before the 
formation of the gaseous disk, and that subsequent spin-up by gas accretion 
from the disk plane has been minimal [134]. 

Coalescence of two black holes during a merger should result in a "spin- 
flip," a reorientation of the spin axis of the more massive black hole. In the 
test-particle limit, the reorientation angle is 

se = H LS oVi-v 2 + o(q 2 ) (5i) 

with fi the cosine of the angle between the orbital angular momentum vector 
and the spin axis of the larger hole. When q > 0.2, the spin orientation is 
overwhelmed by the plunging body in a retrograde merger, even if the initial 
spin of the larger hole was close to maximal. Hence, even "minor mergers" 
(defined, following galactic dynamicists, as mergers with q < 0.3) are able to 
produce a substantial reorientation. In fact there is a class of active galaxies 
which exhibit radio lobes at two, nearly-orthogonal orientations, and in which 
the production of plasma along the fainter lobes appears to have ceased [115, 
32]. These "X-shaped" or "winged" radio galaxies, of which about a dozen 
are known, are plausible sites of recent (within the last ~ 10 8 yr) black-hole 
coalescence [137, 233]. Furthermore the implied coalescence rate is roughly 
consistent with the expected merger rate for the host galaxies of luminous radio 
sources [137]. Alternative models have been proposed for the X-shaped sources, 
including a warping instability of accretion disks [174], backflow of gas along 
the active lobes [116], and binary-disk interactions before coalescence [120]. 5 
It is likely that all of these mechanisms are active at some level and that the 
time scale for realignment influences the radio source morphology, with the 
most rapid realignments producing the classical X-shaped sources, while slower 
realignment would cause the jet to deposit its energy into a large volume, leading 
to an S-shaped FRI radio source [137]. 

If the black holes are spinning prior to coalescence, they will experience spin- 
orbit precession, on a time scale that is intermediate between t gr and the orbital 

5 Liu (2004) criticized the black hole coalescence model on the grounds that "calculations 
based on general relativity show that the change in inclination of a rotating central SMBH 
is negligible in a minor merger and a significant reorientation of the active SMBH requires 
a comparatively rare major merger (Hughes & Blandford 2003)." This erroneous statement 
probably had its origin in the final sentence of the Hughes & Blandford paper, which states 
that "An abrupt change in inclination... requires a comparatively rare major merger." Hughes 
& Blandford defined a "major merger" as having a mass ratio q > 0.1, in conflict with 
the standard definition among galactic dynamicists, q > 0.3. In fact Hughes & Blandford 
conclude, in agreement with Merritt & Ekers (2002), that mass ratios exceeding ~ 0.2 can 
result in spin-flips. 
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period. To PN2.5 order, the spin angular momentum of either hole evolves as 
[98] 



dSi 1 



x Si (52) 



where h is a unit vector in the direction of the displacement vector between 
the two black holes. The evolution equation for S2 is given by interchanging 
the indices. The magnitude of each spin vector remains fixed (to this order), 
and each spin precesses around the total angular momentum vector J = L or b + 
Si + S2. When the two black holes are comparably massive, the orbital angular 
momentum greatly exceeds the spin angular momentum of cither hole until just 
prior to coalescence. As the binary shrinks, the spins have a tendency to unalign 
with L or b- Ignoring spin-spin effects, the precession rate for equal-mass holes in 
a circular orbit is 

n 1/2 

"1 = 2-5575- (53) 
The precession rate is lower than the orbital frequency: 

fl OI b 2 a 16 a 



Oi 1 [i 7 R 
where R = 2mi 2 , but higher than the radiation reaction time scale: 



(54) 



1 35 a 3 ' 2 _ 35V2 /a\3/2 

W htgr - W^ZJT 2 - 512T \RJ ' (55) 



m 12 

unless the binary is close to coalescence. If the dominant source of energy loss 
during the late stages of infall is gravitational radiation, the spin direction will 
undergo many cycles of precession before the black holes coalesce. This does not 
seem to happen in the X-shaped radio sources, based on the apparently sudden 
change in jet direction; if the coalescence model for the X-sources is correct, the 
final stages of infall must occur on a shorter time scale than t gr . This might be 
seen as evidence that shrinkage of the binary is usually driven by gas dynamics, 
not gravitational radiation losses, prior to the final coalescence. 

When the black hole masses are very different, q <C 1, the ratio of spin of 
the larger hole to L OID is 



q V a J q \ a 



Lorb 



(56) 



where Si = fim 2 . The two quantities are approximately equal when the sep- 
aration measured in units of the larger hole's Schwarzschild radius is equal to 
q~ 2 . When this separation is reached, the binary orbit rapidly changes its plane, 
and a new regime is reached where the spin of the smaller black hole precesses 
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about the spin of the larger hole. The precession rate of the larger hole is given 

by 

2m/ m 2 

in both regimes, and 

fii i?' 2tt 9r 128tt li?/ v ; 

The spin direction of a black hole formed via binary coalescence is also 
affected by torques that reorient the binary prior to coalescence. The role of 
torques from gaseous accretion disks was discussed above; another source of 
torques is perturbations from passing stars or gas clouds [134]. A single star 
that passes within a distance <~ 3a of the binary will exchange orbital angular 
momentum with it, leading both to a change in the binary's orbital eccentricity 
as well as a change in the orientation of the binary's spin axis, as discussed 
above. Referring to equations (17) and (20), the reorientation rate is related to 
the hardening rate via 

< A ^> = ^Cden (59) 

where t^l Tden = a(d/dt)(l/a). Scattering experiments [134] give L/H w 4 for a 
hard, equal-mass binary. The implied change in the binary's orientation after 
shrinking from a w ah to a w 10 _3 a^ is 



f*«J*°^. (60) 

V TO12 

The reorientation begins to be significant if m Jr /m 12 > 10~ 3 , which may be the 
case for intermediate-mass black holes. 



10 Summary 

1. No completely convincing case of a bound, binary SBH has yet been disc- 
covered. "Dual" SBHs, i.e. two, widely-separated SBHs in a single system, 
have been seen in some interacting galaxies and binary quasars. However 
none has a projected separation less than ~ 1 kpc. Strong limits can be 
placed on the binarity of the Milky Way SBH. 

2. The evolution of binary SBHs in gas-poor galaxies is dominated by the 
gravitational slingshot ejection of stars that pass near the binary, carry- 
ing away energy and angular momentum. Once the binary has ejected 
all stars on intersecting orbits, continued hardening depends on a refilling 
of the binary's "loss cone." Possible mechanisms for loss-cone refilling in- 
clude star-star gravitational scattering, chaotic orbits in non-axisymmetric 
potentials, and perturbations from additional massive objects. 
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3. Most iV-body simulations of binary evolution have been based on such 
small particle numbers that the binary's loss cone was refilled at a spuri- 
ously high rate by gravitational encounters, Brownian motion of the binary 
and other finite-iV effects. As a result, many results from these studies 
can not usefully be extrapolated to the large- N regime; in particular, they 
do not make useful predictions about the binary hardening rates expected 
in real galaxies. 

4. TV-body studies do appear to make robust predictions about the "mass 
deficit," the mass in stars ejected from a galactic nucleus by an evolving 
binary. Observed mass deficits are of order 1 — 2 times the SBH mass, 
consistent with iV-body predictions. 

5. Studies of the interaction of binary SBHs with gas are in their infancy. Ma- 
jor uncertainties are the amount, distribution and thermodynamic state 
of gas very near the centers of galaxies containing massive binaries. 

6. Binary coalescence can have a large influence on SBH spins; even mass 
ratios as extreme as 10 : 1 can substantially spin up or re-orient a SBH. 
Evidence for spin-flips may have been observed in the so-called X-shaped 
radio sources. 
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